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Abstract 

We discuss recent progress in the development of simulation algorithms that do not rely on any concept of quantum theory but 
are nevertheless capable of reproducing the averages computed from quantum theory through an event-by-event simulation. 
The simulation approach is illustrated by applications to Einstein-Podolsky-Rosen-Bohm experiments with photons. 
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1. Introduction 

Computer simulation is widely regarded as comple- 
mentary to theory and experiment [1]. The standard 
approach is to start from one or more basic equations of 
physics and to employ a numerical algorithm to solve 
these equations. This approach has been highly suc- 
cessful for a wide variety of problems in science and 
engineering. However, there are a number of physics 
problems, very fundamental ones, for which this ap- 
proach fails, simply because there are no basic equa- 
tions to start from. 

Indeed, as is well-known from the early days in the 
development of quantum theory, quantum theory has 
nothing to say about individual events [2,3,4]. Recon- 
ciling the mathematical formalism that does not de- 
scribe individual events with the experimental fact that 
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each observation yields a definite outcome is referred 
to as the quantum measurement paradox and is the 
most fundamental problem in the foundation of quan- 
tum theory [3] . 

In view of the quantum measurement paradox, it is 
unlikely that we can find algorithms that simulate the 
experimental observation of individual events within 
the framework of quantum theory. Of course, we could 
simply use pseudo-random numbers to generate events 
according to the probability distribution that is ob- 
tained by solving the time-independent Schrodinger 
equation. However, the challenge is to find algorithms 
that simulate, event-by-event, the experimental obser- 
vations of, for instance, interference without first solv- 
ing the Schrodinger equation. 

This paper is not about a new interpretation or an 
extension of quantum theory. The proof that there ex- 
ist simulation algorithms that reproduce the results of 
quantum theory has no direct implications on the foun- 
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Fig. 1. (color online) Schematic diagram of an EPRB experiment with photons. 



dations of quantum theory: These algorithms describe 
the process of generating events on a level of detail 
about which quantum theory has nothing to say [3,4]. 
The average properties of the data may be in perfect 
agreement with quantum theory but the algorithms 
that generate such data are outside of the scope of what 
quantum theory can describe. 

In a number of recent papers [5,6,7,8,9,10,11,12,13,14, 
we have demonstrated that locally-connected net- 
works of processing units can simulate event-by-event, 
the single-photon beam splitter and Mach-Zehnder in- 
terferometer experiments, universal quantum compu- 
tation, real Einstein-Podolsky-Rosen-Bohm (EPRB) 
experiments, Wheeler's delayed choice experiment and 
the double-slit experiment with photons. Our work 
suggests that we may have discovered a procedure 
to simulate quantum phenomena using event-based, 
particle-only processes that satisfy Einstein's crite- 
rion of local causality, without first solving a wave 
equation. In this paper, we limit the discussion to 
event-by-event simulations of real EPRB experiments. 



2. EPRB experiments 

In Fig. 1, we show a schematic diagram of an EPRB 
experiment with photons (see also Fig. 2 in [18]). The 
source emits pairs of photons. Each photon of a pair 
propagates to an observation station in which it is ma- 
nipulated and detected. The two stations are separated 
'16 s pl]fcially and temporally [18]. This arrangement pre- 
vents the observation at station 1 (2) to have a causal 
effect on the data registered at station 2 (1) [18]. As the 
photon arrives at station i = 1, 2, it passes through an 
electro-optic modulator that rotates the polarization 
of the photon by an angle depending on the voltage ap- 
plied to the modulator. These voltages are controlled 
by two independent binary random number generators. 
As the photon leaves the polarizer, it generates a sig- 
nal in one of the two detectors. The station's clock as- 
signs a time-tag to each generated signal. Effectively, 
this procedure discretizes time in intervals of a width 
that is determined by the time-tag resolution r [18] . In 
the experiment, the firing of a detector is regarded as 
an event. 



2 



As we wish to demonstrate that it is possible to 
reproduce the results of quantum theory (which im- 
plicitly assumes idealized conditions) for the EPRB 
gedanken experiment by an event-based simulation al- 
gorithm, it would be logically inconsistent to "recover" 
the results of the former by simulating nonideal exper- 
iments. Therefore, we consider ideal experiments only, 
meaning that we assume that detectors operate with 
100% efficiency, clocks remain synchronized forever, 
the "fair sampling" assumption is satisfied [19] , and so 
on. We assume that the two stations are separated spa- 
tially and temporally such that the manipulation and 
observation at station 1 (2) cannot have a causal effect 
on the data registered at station 2 (1). Furthermore, 
to realize the EPRB gedanken experiment on the com- 
puter, we assume that the orientation of each electro- 
optic modulator can be changed at will, at any time. 
Although these conditions are very difficult to satisfy 
in real experiments, they are trivially realized in com- 
puter experiments. 

In the experiment, the firing of a detector is regarded 
as an event. At the nth event, the data recorded on 
a hard disk at station i = 1, 2 consists of x n ,i = ±1, 
specifying which of the two detectors fired, the time tag 
t„ t i indicating the time at which a detector fired, and 
the two-dimensional unit vector a Uii that represents 
the rotation of the polarization by the electro-optic 
modulator. Hence, the set of data collected at station 
i = 1,2 during a run of N events may be written as 



Ti = {x n ,i = ±1, t„,i, a„,i|n = 1,...,N}. 



(1) 



In the (computer) experiment, the data {Ti, T2} may 
be analyzed long after the data has been collected [18]. 
Coincidences are identified by comparing the time dif- 
ferences {i n ,i — tm,2\n, m = 1, . . . , iV} with a time win- 
dow W [18]. Introducing the symbol y' to indicate 
that the sum has to be taken over all events that satisfy 
&i — a n ,j for i = 1,2, for each pair of directions ai and 
a 2 of the electro-optic modulators, the number of co- 
incidences C xy = C xy {a.\,a.2) between detectors D x ,\ 
(x — ±1) at station 1 and detectors D y ^ (y = ±1) at 
station 2 is given by 



J2' s x , XnA s y , Xmi2 e(w-\t„,i-t m M), (2) 



where Q(t) is the Heaviside step function. We empha- 
size that we count all events that, according to the 



same criterion as the one employed in experiment, cor- 
respond to the detection of pairs. The average single- 
particle counts and the two-particle average are defined 
by 



£i(ai,a 2 ) = 

£ 2 (ai,a 2 ) 
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respectively. In Eqs. (3) and (4), the denominator is 
the sum of all coincidences. 

For later use, it is expedient to introduce the function 



S(a, b, c, d) = E(a, c) - E(a, d) + E(b, c) + E(b, d) 
and its maximum 



S m ax= max S(a, b,c,d). 

a.b.c.d 



2.1. Analysis of real experimental data 



(5) 
(6) 



We illustrate the procedure of data analysis and the 
importance of the choice of the time window W by an- 
alyzing a data set (the archives Alice.zip and Bob. zip) 
of an EPRB experiment with photons that is publicly 
available [20]. 

In the real experiment, the number of events de- 
tected at station 1 is unlikely to be the same as the 
number of events detected at station 2. In fact, the 
data sets of Ref. [20] show that station 1 (Alice.zip) 
recorded 388455 events while station 2 (Bob. zip) 
recorded 302271 events. Furthermore, in the real 
EPRB experiment, there may be an unknown shift 
A (assumed to be constant during the experiment) 
between the times t n ,i gathered at station 1 and the 
times tm,2 recorded at station 2. Therefore, there is 
some extra ambiguity in matching the data of station 
1 to the data of station 2. 

A simple data processing procedure that resolves this 
ambiguity consists of two steps [22]. First, we make 
a histogram of the time differences t n ,i — t m ,2 with a 
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Fig. 2. (color online) Smax as a function of the time win- 
dow W, computed from the data sets contained in the 
archives Alice.zip and Bob. zip that can be downloaded from 
Ref. [20]. Bullets (red): Data obtained by using the relative 
time shift A = 4 ns that maximizes the number of coin- 
cidences. Crosses (blue): Raw data (A = 0). Dashed line 
at 2\/2: S max if the system is described by quantum the- 
ory (see Section 3). Dashed line at 2: Smax if the system is 
described by the class of models introduced by Bell [21]. 

small but reasonable resolution (we used 0.5 ns). Then, 
we fix the value of the time-shift A by searching for 
the time difference for which the histogram reaches its 
maximum, that is we maximize the number of coinci- 
dences by a suitable choice of A. For the case at hand, 
we find A = 4 ns. Finally, we compute the coinci- 
dences, the two-particle average, and Smax using the 
expressions given earlier. The average times between 
two detection events is 2.5 ms and 3.3 ms for Alice and 
Bob, respectively. The number of coincidences (with 
double counts removed) is 13975 and 2899 for (A = 4 
ns, W = 2 ns) and (A = , W — 3 ns) respectively. 

In Fig. 2 we present the results for S m ax as a func- 
tion of the time window W . First, it is clear that Smax 
decreases significantly as W increases but it is also 
clear that as W — > 0, Smax is not very sensitive to the 
choice of W [22]. Second, the procedure of maximizing 
the coincidence count by varying A reduces the maxi- 
mum value of Smax from a value 2.89 that considerably 
exceeds the maximum for the quantum system (2\/2, 
see Section 3) to a value 2.73 that violates the Bell in- 
equality (Smax < 2, see Ref. [21]) and is less than the 
maximum for the quantum system. 

Finally, we use the experimental data to show that 
the time delays depend on the orientation of the po- 
larizer. To this end, we select all coincidences between 
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Fig. 3. (color online) Normalized coincidence counts as 
a function of time tag difference t n ,i — t m ,2, computed 
from the data sets contained in the archives Alice.zip and 
Bob. zip [20], using the relative time shift A = 4 ns that 
maximizes the number of coincidences. Bullets (red): 9± = 
and 82 = tt/8; Crosses (blue): 9\ = and 82 = 3tt/8. 



and -D+,2 (see Fig. 1) and make a histogram of 
the coincidence counts as a function of the time-tag dif- 
ference, for fixed orientation 9% = and the two orien- 
tations 82 = 7r/8, 3-7t/8 (other combinations give sim- 
ilar results). The results of this analysis are shown in 
Fig. 3. The maximum of the distribution shifts by ap- 
proximately 1 ns as the polarizer at station 2 is rotated 
by 7r/4, a demonstration that the time-tag data is sen- 
sitive to the orientation of the polarizer at station 2. A 
similar distribution of time-delays (of about the same 
width) was also observed in a much older experimen- 
tal realization of the EPRB experiment [23] . The time 
delays that result from differences in the orientations 
of the polarizers is much larger than the average time 
between detection events, which for the data that we 
analyzed is about 30000 ns. In other words, the loss in 
correlation that we observe as a function of increasing 
W (see Fig. 2) cannot be explained by assuming that 
we calculate correlations using photons that belong to 
different pairs. 

Strictly speaking, we cannot derive the time delay 
from classical electrodynamics: The concept of a pho- 
ton has no place in Maxwell's theory. A more detailed 
understanding of the time delay mechanism requires 
dedicated, single-photon retardation measurements for 
these specific optical elements. 
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2.2. Role of the coincidence window W 

The crucial point is that in any real EPR-type exper- 
iment, it is necessary to have an operational procedure 
to decide if the two detection events correspond to the 
observation of one two-particle system or to the obser- 
vation of two single-particle systems. In standard "hid- 
den variable" treatments of the EPR gedanken exper- 
iment [21], the operational definition of "observation 
of a single two-particle system" is missing. In EPRB- 
type experiments, this decision is taken on the basis of 
coincidence in time [23,25,18]. 

Our analysis of the experimental data shows beyond 
doubt that a model which aims to describe real EPRB 
experiments should include the time window W and 
that the interesting regime is W — > 0, not W — ► oo 
as is assumed in all textbook treatments of the EPRB 
experiment. Indeed, in quantum mechanics textbooks 
it is standard to assume that an EPRB experiment 
measures the correlation [21] 

N 

C X y ^ = ^ ^ ^,Itl,l^S,I»,21 (7) 

n=l 

which we obtain from Eq. (2) by taking the limit 
W —> oo. Although this limit defines a valid theoret- 
ical model, there is no reason why this model should 
have any bearing on the real experiments, in particu- 
lar because experiments pay considerable attention to 
the choice of W. In experiments a lot of effort is made 
to reduce (not increase) W [18,22]. 

3. Quantum theory 

According to the axioms of quantum theory [4] , re- 
peated measurements on the two-spin system described 
by the density matrix p yield statistical estimates for 
the single-spin expectation values 

Si(a) = («7i-a) , S 2 (b) = (a 2 -b), (8) 
and the two-spin expectation value 

S(a,b) = (<ri-a<72-b), (9) 

where Oi — (erf, a\ , erf ) are the Pauli spin- 1/2 matrices 
describing the spin of particle i — 1,2 [4], and a and 
b are unit vectors. We have introduced the tilde to 



distinguish the quantum theoretical results from the 
results obtained from the data sets {Ti, T 2 }. 

The quantum theoretical description of the EPRB 
experiment assumes that the system is represented by 
the singlet state |*) = (|i/>i|V> 2 - |V)i|-ff> 2 ) /V^ of 
two spin- 1/2 particles, where H and V denote the hor- 
izontal and vertical polarization and the subscripts re- 
fer to photon 1 and 2, respectively. For the singlet state 
p = |*)(*|, Ei(a) = E 2 {f3) = and 

E(a,P) = -cos2(a - P). (10) 

4. Simulation model 

A concrete simulation model of the EPRB experi- 
ment sketched in Fig. 1 requires a specification of the 
information carried by the particles, of the algorithm 
that simulates the source and the observation stations, 
and of the procedure to analyze the data. In the fol- 
lowing, we describe a slightly modified version of the 
algorithm proposed in Ref. [9], tailored to the case of 
photon polarization. 

Source and particles: The source emits particles 
that carry a vector S„,i = (cos(£ n + (i— l)-7r/2), sin(£ n + 
(i — l)7r/2), representing the polarization of the pho- 
tons that travel to station i = 1 and station i = 2, 
respectively. Note that S n ,i • S Ui 2 = 0, indicating that 
the two particles have orthogonal polarizations. The 
"polarization state" of a particle is completely charac- 
terized by f n , which is distributed uniformly over the 
whole interval [0, 2ir[. For the purpose of mimicing the 
apparent unpredictability of the experimental data, we 
use uniform random numbers. However, from the de- 
scription of the algorithm, it will be clear that the use 
of random numbers is not essential. Simple counters 
that sample the intervals [0, 27r[ in a systematic, but 
uniform, manner might be employed as well. 

Observation station: The electro-optic modula- 
tor in station i rotates S n ,i by an angle y n ,i, that is 
a„,i = (cos7 nj i,sin7 ni i). The number M of different 
rotation angles is chosen prior to the data collection 
(in the experiment of Weihs et al., M = 2 [18]). We use 
2M random numbers to fill the arrays (ai , (Xm) and 
(Pi, ...,Pm)- During the measurement process we use 
two uniform random numbers 1 < m, m' < M to select 
the rotation angles 7„,i = ct m and 7^2 = p m > ■ The 
electro-optic modulator then rotates S„,j = (cos(£ n + 
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(i- l)7r/2), sin(£ n + (j- l)7r/2) by 7 n ,j, yielding S, M = 
(cos(£„ - j n ,i + l)7r/2),sin(£ n - j n ,i + 1)^/2). 

The polarizer at station i projects the rotated vector 
onto its i-axis: S n ,z ■ Xj = cos(£ n — 7n.i + (i — l) 7r /2), 
where x, denotes the unit vector along the :r-axis of the 
polarizer. For the polarizing beam splitter, we consider 
a simple model: If cos 2 (£ n — j n j + (i — l)ir/2) > 1/2 
the particle causes -D+i,i to fire, otherwise D-i t i fires. 
Thus, the detection of the particles generates the data 
x n ,i = sign(cos2(£„ - y nii + (i - l)n/2)). 

Time-tag model: To assign a time-tag to each 
event, we assume that as a particle passes through 
the detection system, it may experience a time delay. 
In our model, the time delay t n ,i for a particle is as- 
sumed to be distributed uniformly over the interval 
[to, to + T], an assumption that is not in conflict with 
available data [22] . In practice, we use uniform random 
numbers to generate t n ,i- As in the case of the angles 
£„, the random choice of t n ,i is merely convenient, not 
essential. From Eq. (2), it follows that only differences 
of time delays matter. Hence, we may put to = 0. The 
time-tag for the event n is then t n ,i £ [0, T]. 

There are not many options to make a reasonable 
choice for T. Assuming that the particle "knows" its 
own direction and that of the polarizer only, we can 
construct one number that depends on the relative an- 
gle: S n ,i-Xi. Thus, T = T(£ n — Jn,i) depends on £„— 7„,i 
only. Furthermore, consistency with classical electro- 
dynamics requires that functions that depend on the 
polarization have period n [24]. Thus, we must have 
T(£n - j n ,i + (i- 1W2) = F((S n ,< ■ xO 2 ). We already 
used cos 2(£ n — y n ,% + (i — l) 7r /2) to determine whether 
the particle generates a + 1 or — 1 signal. By trial and er- 
ror, we found that T(f n -0i) = T F(\ sin2(£ n -0i)|) = 
To | sin2(£„ - 6»i)| d yields useful results [9,10,11,12,13]. 
Here, To = maxe T(6) is the maximum time delay and 
defines the unit of time, used in the simulation and d is 
a free parameter of the model. In our numerical work, 
we set To — 1. 

Data analysis: For fixed iV and M, the algorithm 
generates the data sets Ti just as experiment does [18]. 
In order to count the coincidences, we choose a time- 
tag resolution < r < To and a coincidence window 
t < W . We set the correlation counts C xy (a m , P m ') to 
zero for all x,y = ±1 and m,m' — 1, M. We com- 
pute the discretized time tags k n ,i = \t n ,i/T~\ for all 
events in both data sets. Here \x~\ denotes the small- 
est integer that is larger or equal to x, that is \x~\ — 
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Fig. 4. (color online) Comparison between computer simu- 
lation data (red bullets) and quantum theory (black solid 
line) for the two-particle correlation E(a,/3). 

1 < a; < \x]. According to the procedure adopted 
in the experiment [18], an entangled photon pair is 
observed if and only if |fe n ,i — fcn,a| < k — \W/t~\. 
Thus, if \k n ,l — knfil < k, we increment the count 

C;c n) 1,^:71,2 {phn i Pm / ) • 

5. Simulation results 

The simulation proceeds in the same way as the ex- 
periment, that is we first collect the data sets { Yi , Y2}, 
and then compute the coincidences Eq. (2) and the cor- 
relation Eq. (4). The simulation results for the coinci- 
dences C xy (a, (5) depend on the time-tag resolution r, 
the time window W and the number of events N, just 
as in real experiments [18] . 

Figure 4 shows simulation data for E(a, f3) as ob- 
tained for d=2, N = 10 6 and W = t = 0.00025T . In 
the simulation, for each event, the random numbers 1 < 
A n: i < M select one pair out of {(on, /3j)\i, j = 1,M}, 
where the angles on and f3j are fixed before the data 
is recorded. The data shown has been obtained by al- 
lowing for M = 20 different angles per station. Hence, 
forty random numbers from the interval [0,360[ were 
used to fill the arrays (on, . . . , otu) and (/3i, . . . , /3m). 
For each of the N events, two different random num- 
ber generators were used to select the angles ot m and 
P m i . The statistical correlation between m and m' was 
measured to be less than 10~ 6 . 

From Fig. 4, it is clear that the simulation data for 
E(a, P) are in excellent agreement with quantum the- 
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ory. Within the statistical noise, the simulation data 
(not shown) for the single-spin expectation values also 
reproduce the results of quantum theory. 

Additional simulation results (not shown) demon- 
strate that the kind of models described earlier are ca- 
pable of reproducing all the results of quantum the- 
ory for a system of two S=l/2 particles [9,10,11,12,13]. 
Furthermore, to first order in W and in the limit that 
the number of events goes to infinity, one can prove 
rigorously that these simulation models give the same 
expressions for the single- and two-particle averages as 
those obtained from quantum theory [9,10,11,12,13]. 



6. Discussion 

Starting from the factual observation that experi- 
mental realizations of the EPRB experiment produce 
the data {Ti, Y2} (see Eq. (1)) and that coincidence 
in time is a key ingredient for the data analysis, we 
have described a computer simulation model that sat- 
isfies Einstein's criterion of local causality and, exactly 
reproduces the correlation i?(ai,a2) = — ai • a.2 that 
is characteristic for a quantum system in the singlet 
state. 

We have shown that whether or not these simulation 
models produce quantum correlations depends on the 
data analysis procedure that is performed (long) after 
the data has been collected: In order to observe the 
correlations of the singlet state, the resolution r of the 
devices that generate the time-tags and the time win- 
dow W should be made as small as possible. Disregard- 
ing the time-tag data (d — or W > To) yields results 
that disagree with quantum theory but agree with the 
models considered by Bell [21]. Our analysis of real ex- 
perimental data and our simulation results show that 
increasing the time window changes the nature of the 
two-particle correlations [9,10,11,12,13]. 

According to the folklore about Bell's theorem, a 
procedure such as the one that we described should 
not exist. Bell's theorem states that any local, hidden 
variable model will produce results that are in conflict 
with the quantum theory of a system of two S = 1/2 
particles [21]. However, it is often overlooked that this 
statement can be proven for a (very) restricted class 
of probabilistic models only. In fact, Bell's theorem 
does not necessarily apply to the systems that we are 



interested in as both simulation algorithms and actual 
data do not need to satisfy the (hidden) conditions 
under which Bell's theorem hold [26,27,28,29]. 

Furthermore, the apparent conflict between the fact 
that there exist event-based simulation models that 
satisfy Einstein's criterion of local causality and repro- 
duce all the results of the quantum theory of a system 
of two S = 1/2 particles and the folklore about Bell's 
theorem, stating that such models are not supposed 
to exist dissolves immediately if one recognizes that 
Bell's extension of Einstein's concept of locality to the 
domain of probabilistic theories relies on the hidden, 
fundamental assumption that the absence of a causal 
influence implies logical independence [30,31]. 

The simulation model that is described in this pa- 
per is an example of a purely ontological model that 
reproduces quantum phenomena without first solving 
the quantum problem. The salient features of our sim- 
ulation models [5,6,7,8,9,10,11,14] are that they 

(i) generate, event-by-event, the same type of data 
as recorded in experiment, 

(ii) analyze data according to the procedure used in 
experiment, 

(iii) satisfy Einstein's criterion of local causality, 

(iv) do not rely on any concept of quantum theory or 
probability theory, 

(v) reproduce the averages that we compute from 
quantum theory. 
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